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The  recovery  of  energy  from  exhaust  gases  using  thermoelectric  generators  is  of  growing  interest.  The 
electrical  loading  of  the  thermoelectric  system  impacts  the  amount  of  energy  that  can  be  recovered.  This 
paper  presents  a  model  for  the  theoretical  limit  of  electrical  power  generation  that  provides  optimal  elec¬ 
trical  loading  conditions  for  a  given  exhaust  stream  and  system  configuration.  The  analysis  of  a  sample 
heat  recovery  configuration  indicates  that  the  simple  isothermal  modeling  approach  often  applied  to 
individual  thermoelectric  leg  pairs  is  not  sufficient  to  optimize  power  generation  when  a  significant 
amount  of  energy  is  removed  from  the  exhaust  stream  via  a  sequence  of  leg  pairs.  The  analysis  estab¬ 
lishes  that  ZT,  the  thermoelectric  figure  of  merit,  is  not  a  sufficient  metric  to  describe  system  level  per¬ 
formance.  The  analysis  predicts  that  there  is  an  optimum  number  of  thermoelectric  leg  pairs  that 
maximize  the  power  extracted  for  any  system,  and  that  adding  more  leg  pairs  beyond  this  optimum 
can  degrade  system  performance.  The  theoretical  limit  for  power  generation  is  compared  to  proposed 
electrically  loading  strategies  found  in  the  literature.  The  developed  model  has  low  computational  load 
and  is  suitable  for  use  in  system  optimization  models. 

©  2014  Elsevier  Ltd.  All  rights  reserved. 


1.  Introduction 

The  growing  need  for  electricity,  and  the  heightened  concern 
about  environmental  issues  associated  with  traditional  energy 
generating  technologies,  has  increased  interest  in  the  use  of 
thermoelectrics  to  directly  convert  thermal  energy  into  electrical 
energy.  Thermoelectrics  may  enable  a  viable  pathway  for  capturing 
waste  heat,  and  may  be  integrated  with  other  power  cycles 
to  improve  overall  efficiency  and  reduce  costs.  There  is  substantial 
interest  in  utilizing  thermoelectric  devices  for  automobile  exhaust- 
to-electricity  conversion  [1].  According  to  a  recent  Department  of 
Energy  Industrial  Technology  Program  study  [2],  the  waste  energy 
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equivalent  to  1.72  billion  barrels  of  oil  is  lost  each  year  in  the 
United  States  industrial  sector  alone.  Significant  energy  savings 
may  be  achieved  by  recovering  even  a  small  portion  of  this 
discarded  thermal  energy.  Thermoelectrics  have  some  distinct 
advantages  since  they  are  modular,  have  no  moving  parts,  and 
are  able  to  operate  over  a  wide  range  of  transient  temperature  con¬ 
ditions.  High  costs  and  poor  efficiencies,  however,  have  deterred  a 
more  rapid  adoption  of  thermoelectric  devices. 

There  have  been  significant  advancements  in  thermoelectric 
material  properties  in  recent  years  [3,4].  Optimization  of  a  single 
leg  pair  to  extract  the  most  waste  energy  has  been  well  established 
by  altering  material  properties  and  leg  geometry  [5].  The  use  of 
functionally  graded  materials  or  the  stacking  of  leg  materials  has 
also  been  considered  to  optimize  individual  devices  [6,7].  How¬ 
ever,  the  optimal  conditions  at  which  to  operate  a  single  leg  pair 


R.J.  Stevens  et  al./ Applied  Energy  133  (2014)  80-88 


81 


Nomenclature 

A 

cross  sectional  area  of  a  leg  (m2) 

T ambient 

sink  temperature  (K) 

b 

length  of  a  leg  pair  including  open  space  between  legs 

V 

voltage  across  leg  pair  (V) 

(m) 

w 

width  of  a  thermoelectric  leg  pair  couple  including  open 

C 

coefficients  used  to  simplify  expressions 

space  between  legs  (m) 

cP 

specific  heat  (J/kg  K) 

X 

distance  from  leading  inlet  of  system  (m) 

D 

coefficients  used  to  simplify  expressions 

z 

thermoelectric  figure  of  merit  (1/K) 

I 

current  (A) 

j 

current  flux  (A/m2) 

Greek  symbols 

k 

thermal  conductivity  (W/m  K) 

(X 

Seebeck  coefficient  (V/K) 

I< 

thermal  conductance  of  a  leg  pair  (W/K) 

n 

efficiency  of  a  thermoelectric  leg  pair 

l 

height  of  a  leg  (m) 

EJsys 

system  efficiency 

L 

system  length  (m) 

X 

Lagrange  multipliers 

m 

ratio  of  electrical  load  to  internal  resistance  for  local 

p 

electrical  resistivity  (Q  cm) 

efficiency  optimization  () 

function  to  be  optimized 

m 

mass  flow  rate  (g/s) 

N 

number  of  leg  pairs  () 

Siihcrrints 

P 

system  electric  power  generation  (W) 

1-6 

constants  index 

Q 

heat  rate  of  a  leg  pair  (W) 

c 

thermal  contact 

Q 

system  heat  transfer  rate  (W) 

h,c 

hot  and  cold  side 

R 

overall  thermal  resistances  between  the  fluids  and  ther¬ 

i,o 

inlet  and  outlet  of  system 

moelectric  junction  for  a  single  leg  pair  (K/W) 

p,n 

p  and  n  type  legs 

K 

thermal  contact  resistance  (K  m2/W) 

r  j  t'  o 

Re 

electrical  resistance  of  a  leg  pair  (Q) 

T 

temperature  (K) 

superscripts 

//  _ _ 1 _ • 

T 

average  of  hot  and  cold  side  leg  temperature  (K) 

on  died  udsis 

J* 

integral  average  voltage  difference  across  leg  pairs  for 

entire  system  (K) 

does  not  necessarily  apply  to  a  sequence  of  leg  pairs,  especially 
when  the  temperature  of  the  heat  source  and/or  sink  changes 
appreciably  along  the  system.  This  is  not  surprising,  since  energy 
removal  by  a  single  leg  pair  affects  the  temperature  field  for  subse¬ 
quent  leg  pairs,  and  the  interaction  between  these  leg  pairs  is  thus 
nontrivial. 

Miller  et  al.  modeled  a  thermoelectric  system  with  an  organic 
Rankine  bottoming  cycle  for  heat  recovery  from  exhaust  heat  [8]. 
The  model  accounted  for  the  effects  of  various  design  parameters 
and  the  thermal  resistance  of  the  heat  exchanger  surfaces.  The 
model  assumes  an  isothermal  surface,  which  is  applicable  only  if 
a  small  fraction  of  the  available  energy  is  extracted  from  the 
exhaust  heat.  The  isothermal  assumption  has  been  used  in  model¬ 
ing  waste  heat  recovery  and  solar  thermal  cycles  [9,10]. 

Gao  et  al.  developed  a  model  of  a  sequence  of  thermoelectric 
modules  in  a  counter  flow  to  recover  waste  energy  from  a  PEM  fuel 
cell  [11].  The  model  enabled  each  module,  a  discrete  device  with  a 
multiple  leg  pairs  wired  electrically  in  series,  to  either  operate 
electrically  in  series  or  in  stand-alone  mode.  For  the  stand-alone 
mode,  the  applied  electrical  load  across  a  module  was  set  to  equal 
the  internal  electrical  resistance  of  the  module.  Although  this 
approach  will  naturally  guarantee  that  the  power  will  be  maxi¬ 
mized  for  a  given  temperature  difference  across  a  single  leg  pair, 
the  approach  will  not  necessarily  generate  the  maximum  power 
for  the  entire  system.  Zhou  et  al.  implemented  a  multi-scale  ther¬ 
moelectric  model  that  also  matched  the  external  load  with  the 
internal  electrical  resistance,  and  explored  the  impact  of  a  range 
of  design  parameters  and  flow  conditions  on  pressure  drop  and 
power  generation  [12].  Yu  and  Zhao  and  Kumar  et  al.  also  consid¬ 
ered  legs  in  series  and  the  load  resistance  was  selected  to  equal 
that  of  the  effective  internal  resistance  of  the  modules  [13-15]. 

Crane  developed  a  numerical  model  for  multiple  leg  pairs  and 
used  a  MATLAB  function  to  iteratively  optimize  each  element  effi¬ 
ciency  for  an  automobile  exhaust  recovery  system  [16].  The  legs 
operated  electrically  in  series.  Crane’s  approach  adjusted  the  total 


load  resistance  in  order  to  optimize  overall  power  generation  as 
opposed  to  the  approach  of  load  matching  typically  done  with  a 
single  module  (a  discrete  device  with  multiple  leg  pairs  in  series). 
He  validated  the  model  experimentally  for  a  small  system.  Wang 
et  al.  used  a  numerical  approach  to  optimize  a  thermoelectric  gen¬ 
erator  coupled  to  vehicle  exhaust,  and  found  that  the  optimal 
power  is  reached  when  the  applied  electrical  resistance  is  greater 
than  the  internal  resistance  for  cases  where  all  leg  pairs  are  in  ser¬ 
ies  [17]. 

Liang  et  al.  used  a  slightly  different  approach  to  optimize  a 
sequence  of  thermoelectric  modules  by  fixing  the  voltage  across 
parallel  strings  of  thermoelectric  legs  pairs  [18].  Min  and  Rowe 
[19]  assumed  that  the  optimal  operation  occurs  when  all  leg  pairs 
individually  operate  at  peak  efficiency  as  defined  by  Ioffe  [20]. 

The  above-cited  literature  indicates  that  the  electrical  loading  of 
a  sequence  of  thermoelectric  elements  has  an  impact  on  system 
overall  performance.  The  purpose  of  this  paper  is  to  develop  a  the¬ 
oretical  limit  for  the  amount  of  power  that  can  be  recovered  from 
an  exhaust  gas.  This  limit  is  extracted  from  a  model  of  low  compu¬ 
tational  load,  so  it  is  suitable  for  use  to  optimize  large  systems.  The 
results  of  this  theoretical  limit  model  are  compared  to  alternative 


Fig.  1.  Schematic  of  thermoelectric  power  generation  system. 
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approaches  to  optimization,  some  of  which  have  been  cited  above. 
In  particular,  we  compare  our  approach  and  results  for  the  follow¬ 
ing  cases:  (i)  the  local  efficiency  of  each  leg  pair  is  optimized,  (ii)  a 
constant  voltage  is  applied  across  all  leg  pairs,  and  (iii)  a  constant 
current  is  applied  across  all  leg  pairs. 

2.  Thermal  model 

A  schematic  of  the  waste  heat  thermoelectric  power  generator 
to  be  modeled  is  shown  in  Fig.  I.  This  structure  is  typical  of  many 
waste  heat  recovery  systems.  A  hot  fluid  flows  across  a  surface, 
which  is  in  thermal  contact  with  multiple  thermoelectric  leg  pairs. 
The  underside  of  these  thermoelectric  leg  pairs  are  cooled  by  cool¬ 
ant  passing  across  a  surface  in  thermal  contact  with  the  bottom 
side  of  the  thermoelectric  leg  pairs.  One  or  both  surfaces  often 
include  a  system  of  fins  to  enhance  heat  transfer.  A  portion  of 
the  energy  transferred  to  the  hot  side  as  heat  is  directly  converted 
to  electrical  energy,  while  the  remaining  energy  is  dumped  to  the 
coolant.  As  energy  is  extracted  from  the  hot  fluid  and  dumped  to 
the  cold  fluid,  the  temperature  distribution  along  the  length  of 
the  system  in  both  the  hot  and  cold  fluids  will  change. 

The  temperatures  of  the  hot  and  cold  fluids  are  defined  as  Th 
and  Tc  respectively  and  vary  with  respect  to  x,  the  distance  along 
the  length  of  the  system;  the  subscripts  h  and  c  will  be  used  gen¬ 
erally  in  what  follows  to  refer  to  the  hot  and  cold  fluids,  respec¬ 
tively.  The  mass  flow  rates,  m,  are  defined  as  being  positive  for 
both  the  hot  and  cold  fluids  moving  in  the  positive  x-direction 
(i.e.,  co-current  flow  as  shown  in  Fig.  1 ).  The  countercurrent  config¬ 
uration  is  modeled  by  choosing  opposite  signs  in  the  mass  flow 
rates.  The  total  rate  of  energy  transferred  to  all  thermoelectric  ele¬ 
ments  in  the  form  of  heat  is 

Qh  —  (tflCp)h(T/!  —  Tho)  (1) 

where  cp  is  the  heat  capacity  of  the  fluid  and  the  subscripts  i  and  o 
denote  the  system  inlet  and  outlet  temperatures,  respectively.  The 
total  rate  of  energy  transferred  from  the  thermoelectric  elements 
on  the  coolant  side  of  the  system  is 

Q-c  —  (iticp)c(Tc,0  —  Tcj)  (2) 

The  electrical  power  generated  by  the  entire  system  is  the  dif¬ 
ference  in  the  two  energy  rates  in  Eqs.  (1)  and  (2) 

P  =  dk-dc  (3) 

One  goal  in  designing  a  heat  recovery  system  is  often  to 
maximize  the  electrical  power  generated,  P,  compared  with  the 
maximum  rate  of  thermal  energy  available  in  the  waste  stream, 
Qmax •  The  system  efficiency  quantifies  the  ratio  of  the  two  as 

llsys  =  P  /  Qmax  (4) 

where 

Qmax  —  (WCp)h(Th,i  ~  T ambient )  (4b) 

In  Eq.  (5),  Tambient  is  the  temperature  of  the  local  heat  sink  which 
is  assumed  in  this  paper  to  be  that  of  the  coolant  fluid  at  its  inlet  to 
the  system,  Tc>i.  This  definition  of  system  efficiency  was  also  used 
by  Zhou  et  al.  [12].  Note  this  definition  of  efficiency  is  different 
from  the  usual  definition  of  efficiency,  p,  defined  as  the  ratio  of 
the  power  generated  by  all  leg  pairs  to  the  heat  absorbed  at  the 
hot  side  of  all  leg  pairs  and  expressed  as 

V  =  P/Qh  (5) 

Eq.  (5)  is  appropriate  for  closed  system  heat  sources  such  as  radio¬ 
isotope  thermoelectric  generators.  Most  waste  heat  recovery  appli¬ 
cations  are  open  systems  so  it  is  more  appropriate  to  consider  the 
ratio  of  generated  power  to  thermal  power  supplied  as  discussed 


by  Min  and  Rowe  [19].  Eq.  (4)  is  slightly  different  than  defined  by 
Min  and  Rowe,  however;  in  Eq.  (4)  the  denominator  is  the  thermal 
power  available  versus  thermal  power  supplied.  The  Min  and  Rowe 
definition  is  more  appropriate  when  fuel  is  being  supplied  solely  for 
the  purposes  of  thermoelectric  power  conversion.  For  applications 
where  waste  heat  is  being  tapped,  Eq.  (4)  is  most  appropriate. 

To  determine  the  temperature  profiles  of  the  fluids  and  junc¬ 
tion,  the  thermal  transfer  at  the  thermoelectric  elements  must  be 
examined.  Consider  a  single  leg  pair  as  shown  in  Fig.  2,  where 
the  n-type  leg  is  immediately  behind  p- type  leg.  The  rate  of  heat 
transfer  to  the  hot  side  of  the  leg  pair  can  be  expressed  as  [5] 

qh  =  ahI  +  K(h-T2)-^  (6) 

where  junction  temperatures  for  the  hot  and  cold  side  of  the 
thermoelectric  elements  are  Ti  and  I2,  respectively  and  also  vary  with 
respect  to  x  (Fig.  2 ).  The  electrical  resistance,  Re,  thermal  conductance, 
K,  and  Seebeck  coefficient,  a,  for  the  leg  pair  can  be  expressed  as 

n  PJn  Pplp  1  K  ?  K 

R'=l^  +  J^  +  2An+2Ap 

rr  _  kpAp  (7) 

A  ~  I  H  / 

Ln  lp 

oc  =  ocp  —  ocn 

where  l  is  the  leg  lengths,  A  is  the  cross  sectional  areas  of  the  legs,  R" 
is  the  electrical  contact  resistance,  p  is  the  electrical  resistivity  of 
the  legs,  and  k  is  the  thermal  conductivity,  all  of  which  are  assumed 
to  be  constant.  In  Eq.  (7)  and  throughout  this  paper,  the  subscripts  n 
and  p  on  any  parameter  denote  an  n-type  or  p- type  semiconductor, 
respectively.  It  is  assumed  that  the  process  is  at  steady  state,  and 
the  surfaces  of  each  thermoelectric  leg  are  isothermal  and  that  they 
have  constant  properties.  The  Thomson  effect  can  be  accounted  for 
by  using  an  integral  average  Seebeck  coefficient  as  demonstrated  by 
Sandoz-Rosado  et  al.  [21]. 

A  similar  expression  for  the  cold  side  of  the  junction  can  be 
expressed  as  5] 

qc  =  &T2I  +  K(T\  - 12)  +~y~  (8) 


The  heat  transfer  to  and  from  the  two  sides  of  the  leg  pair  can 
also  be  expressed  in  terms  of  the  temperature  differences  between 
the  source  and  sink  fluids  and  thermal  resistances  between  the 
fluids  and  thermoelectric  junctions  as 


(9a) 

(9b) 


Tc  (x)  •  Cold  Fluid 

- b  — 


l 


I'M 

R, 


Rc. 

r.(*) 

Re 


Fig.  2.  Schematic  of  thermoelectric  leg  with  temperature  definitions. 
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where  Rh  and  Rc  are  the  overall  thermal  resistances  between  the  flu¬ 
ids  and  thermoelectric  junction  for  a  single  leg  pair.  These  resis¬ 
tances  account  for  the  ceramic  layer  typically  used  in  modules  as 
well  as  the  fins  systems  and  convection  resistances  (Fig.  2). 

Although  there  are  discrete  modules  that  comprise  the  thermo¬ 
electric  system  in  Fig.  1,  we  will  model  the  system  in  such  a  way 
that  the  system  can  be  viewed  as  a  single  continuous  thermoelec¬ 
tric  from  which  a  current  can  be  drawn  locally.  To  this  end,  we 
identify  a  length  and  width  of  the  thermoelectric  leg  pair,  denoted 
respectively  as  b  and  w,  that  extend  beyond  the  physical  length  and 
width  of  the  pair,  adding  in  half  of  the  spaces  to  the  adjacent  leg 
pairs  shown  in  Fig.  2. 

We  then  assume  that  the  length  of  any  leg  pair  is  small  com¬ 
pared  with  the  total  length  of  the  system,  i.e.,  h/L<  1,  which 
allows  us  to  view  the  system  as  a  single  continuous  module  along 
which  the  temperatures  and  current  can  vary.  Consistent  with  this 
approximation,  we  define  a  local  current  flux  (i.e.,  current  per 
area),  j(x).  The  current  in  the  kth  module  can  be  expressed  as: 

rxk+ 1 

h=w  j(x)dx  (10a) 

Jxk 

where 


Cl  =  w/(mcp)hR l 
C2  =  w/(mcp)cR" 
Di  =  aR|| 

D2  =  I<"K 
D3=R'X/ 2 
D4  =  aR" 

D5  =  K"R" 

D6  =  R'X/2 


(13e) 


In  (13e),  note  that  the  mass  flow  rates,  rh,  are  for  a  single  row  of 
thermoelectric  leg  pairs  with  a  width  of  w,  and  D2  and  D5  are  the 
ratios  of  the  fluid-thermoelectric  junction  resistance  to  the  thermal 
resistance  across  the  thermoelectric  elements.  Eq.  (13)  are  subject 
to  the  boundary  condition  Th(x  =  0)  =  Thti ,  which  is  the  temperature 
of  the  incoming  waste  heat  fluid.  The  boundary  condition  for  the 
cold  side  is  either  Tc(x  =  0)  =  TCti  for  co-current  flow  or  Tc(x  =  L)  =  TCfi 
for  countercurrent  flow. 

Eq.  (13)  define  the  temperature  profiles  for  an  arbitrary  current 
flux  profile,  j(x).  The  power  that  can  be  extracted  is  obtained  by 
subtracting  Eq.  (lib)  from  Eq.  (11a)  and  integrating  over  the  entire 
length  of  the  system  to  obtain 


xfc  =  (k-l)b,ke[l,N  +  l]  (10b) 

In  (10b),  N  is  the  number  of  leg  pairs.  We  note  that  in  the  limit  as 
b/L  <  1,  we  can  approximate  (10)  as: 

/~wbj(x)  (10c) 

where  the  current  is  now  interpreted  as  being  continuous. 

In  accordance  with  definition  (10b)  and  the  geometry  shown  in 
Fig.  2,  L  =  b-N.  To  make  use  of  Eq.  (10)  all  parameters  and  variables 
are  placed  on  a  per  area  basis,  and  to  this  end,  Eqs.  (6),  (8),  (9)  are 
written  as: 


p  =  w  [  [a(T,  -  T2)j  -  R'f]  dx  (14) 

Of  interest  is  to  maximize  the  total  power  generated  in  accor¬ 
dance  with  Eq.  (14)  so  that  the  system  efficiency  given  by  Eq.  (4) 
is  maximized;  this  is  obtained  by  determining  an  optimal  current 
flux  profile,  j(x).  The  determination  of  this  current  flux  is  described 
in  the  next  section. 

3.  Power  optimization— the  theoretical  limit  model 


(Th  =  ~ ^  =  «hj  +  K"(h 

(31a) 

T  T  R"  P 

<?c  =  —ft”  ~  +  K”(h 

(lib) 

p"  =  c/JT]  -  T2)j  -  Rf 

(11c) 

where 

K  =  Ri,wb 

R!'c  =  Rcwb 
fi"  =  Rewb 

K"  =  K/wb 

(lid) 

The  heat  extracted  from  the  hot  fluid  and  rejected  to  the  cold 


fluid  can  be  expressed  as 

Th\  —  Th  L 

q'h  =  -{mcp)h — -f1 - t~ 

p,h  w(xi+1  -  Xj) 

(12a) 

TCL  -  Tc  L 

(12b) 

Combining  Eqs.  (11)  and  (12)  while  considering  the  limit  of 
b/L  <c  1,  the  following  system  of  equations  is  obtained 


^  =-C,(T„-Tt)  (13a) 

^  =  C2(T2-Tc)  (13b) 

Th  -  T,  =  D,  r,j  +  D2  (7,  -  T2)  -  D3f  (13c) 

T2-Tc  =  DtTzj  +  D5  (T,  -T2)+  Dej2  (13d) 

where 


The  appearance  of  the  current  density  in  the  integral  expression 
for  the  power,  Eq.  (14),  suggests  that  a  variational  method  may  be 
used  to  maximize  the  power.  To  this  end,  we  rewrite  Eqs.  (13c,d)  to 
obtain  two  explicit  expressions  for  Ti(x)  and  T2(x)  in  terms  of  Th(x), 
Tc(x)  and  j(x)  as 


j  .  x  _D2Tc  +  (1  +D5)Tfr  -D4]Th  +  (D2D6+D3D5+D3f-D3D4j3 
(1+D2+D5)  +  (D1D5-D2D4  +  D1-D4)j-D1D4/2 
j  ,  '  D5T1,  +  (1  +D2)7c  +  DijTc  +  (D2D6  +  D3D5  +  Def_+D^D6f  (151-,) 

2 (1+D2+D5)  +  (D1D5-D2D4  +  D1-D4y-D1D4;-2 

The  total  power  in  Eq.  (14)  is  thus  seen  to  be  a  function  of  Th(x ), 
Tc(x)  and  j(x).  We  need  to  maximize  the  power,  P  in  Eq.  (14),  sub¬ 
ject  to  constraints  imposed  by  Eqs.  (13a)  and  (13b).  To  do  so,  a 
Lagrangian  function  i Jj  is  constructed  based  on  the  integrand  from 
the  power  expression  (14),  Eqs.  (13a)  and  (13b),  and  x-dependent 
Lagrange  multipliers  2i(x)  and  22(x)  to  yield:  [22] 


=  [«(Ti  -  T2)j  -  R'f]  +  h  (x)  [T'h  +  C,  (Th  -  T, )] 

+  A2(x)[T'-C2(T2-Tc)]  (16) 

where  T'h  =  ^  and  T'c  =  By  inspection,  we  see  that  when  the 
constraints  given  in  Eqs.  (13a)  and  (13b)  are  satisfied,  the  terms 
multiplying  the  Lagrange  multipliers  are  zero,  and  the  power  is 
given  as 


P  = 


(17) 


Optimization  of  the  functional,  i/q  with  respect  to  the  three 
functions  Th(x),  Tc(x)  and  j(x),  yields  the  well-known  Euler- 
Lagrange  equations  [22];  system  specific  boundary  conditions  are 
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necessary  to  eliminate  boundary  terms  arising  through  the  stan¬ 
dard  integration  by  parts  in  the  derivation  of  the  equations  [22]. 
These  equations  and  constraints  are  given  as: 

i-° 


3±_±(3t\=0 

dTh  dx  \OT'J 
d< 1/  d  (d\\> 


3TC  dx  \dTc )  ° 

@x  =  0,  =  0 

d  T'c 

@x  =  L,  =  0 

dT'h 

The  substitution  of  Eq.  (16)  into  Eq.  (18)  yields  the  following 

(dTx  dT2 


(18) 


dj  dj 
dT1_dT1 
dTh  dTh 
dT1_dT1 
dTc  dTc 


+  a(T,  -  T2)  -  2jR"  -  C,A,  ^  -  C2A2  ^  =  0 


+  fa 
+  fa 


C  -C 
1  Cl  9Th_ 

C2~C2W 

dTc 


i)T2  dAi  _  n 

-C2X2Wh-^~° 

dT-t  d/,2 

-C^Wc~W  =  0 


(19a) 


with  the  following  boundary  conditions. 


@x  =  0,  Th  =  Thm.  A2  =  0  ;  @x  =  L,  Tc  =  Tc4n,  X\  =  0. 

(19b) 

The  set  of  nonlinear  Eq.  (19),  coupled  with  Eqs.  (13)  and  (15), 
represent  a  system  of  5  equations  in  the  5  unknown  functions: 
Th(x),  Tc(x),  j(x),  fa(x)  and  fa(x);  this  system  is  solved  using  a  full 
Newton’s  iteration  method.  Sufficient  accuracy  is  obtained  with 
400  equally  spaced  discretizations  along  x  for  each  variable.  The 
coupled  system  is  solved  iteratively  until  both  the  Euclidean  norm 
of  the  equation  residuals  and  that  of  the  change  to  all  dependent 
variables  is  less  than  10-8. 

In  the  case  of  co-current  flow,  the  boundary  conditions  in  the 
system  (19)  switch  location  in  such  a  way  that  there  are  two 
boundary  conditions  at  each  end  of  the  computational  domain. 
The  boundary  conditions  for  co-current  flow  are 

@  X  =  0,  Th  =  Thin,  Tc  =  Tcin;  and  @  x  =  L,  fa  =  0,  fa  =  0 

(19c) 

For  purposes  of  discussion,  we  will  refer  to  the  above  described 
approach  as  the  Theoretical  Limit  Model. 


4.  Results  and  discussion 

We  now  examine  results  of  the  theoretical  limit  model.  To  do 
so,  we  specify  some  typical  properties  (Table  1 )  and  system  level 
parameters  (Table  2)  that  will  be  used  for  all  examples  in  the  rest 
of  this  paper,  except  where  noted.  With  this  initial  specification, 
we  consider  a  waste  heat  recovery  application  with  an  exhaust 
temperature  of  500  °C  and  mass  flow  rate  of  25  g/s  (these  values 
indicated  in  Table  2).  The  system  is  to  be  cooled  in  a  counter  flow 
configuration  by  a  coolant  at  50  °C  with  a  heat  capacity  rate,  mcp, 
ten  times  that  of  the  exhaust  side. 

In  Table  2,  we  have  adhered  to  the  naming  convention  that  a 
module  corresponds  to  multiple  leg  pairs  (modules  are  how  ther¬ 
moelectric  are  used  in  practice);  however,  note  that  it  is  the  perfor¬ 
mance  of  each  leg  pair  in  the  system  that  is  being  examined  in  our 
work.  The  model  developed  in  Sections  2  and  3  is  used  to  extract 
the  optimal  current  flux  profile,  j(x),  that  should  be  applied  to 


Table  1 

Leg  pair  parameters. 


Parameters 

P  type 

n  type 

Seebeck  coefficient,  a 

200  10-6  V/K 

— 15010-6  V/K 

Electrical  resistivity,  p 

1.510-3  Q  cm 

2.010-3  Q  cm 

Thermal  conductivity,  k 

1.1  W/m  K 

1.3  W/m  K 

Leg  length,  l 

1  mm 

1  mm 

Leg  area,  A 

2.25  mm2 

2.25  mm2 

Electrical  contact  resistance,  R" 

0  Q  m2 

OQm2 

Table  2 

Other  system  parameters. 


Exhaust  inlet  temperature,  Thiin  500  °C 

Coolant  inlet  temperature,  TCjin  50  °C 

Exhaust  mass  flow  rate,  m  25  g/s 

Coolant  heat  capacity  ratio  to  exhaust  heat  capacity  ratio,  1 0 

(mcp)c/(mcp)h 

#  Leg  pairs  per  module  127 

#  Modules  in  system  25 

Ratio  of  hot  side  external  resistance  to  internal  thermal  0.1 

resistance,  I<"R ^ 

Ratio  of  cold  side  external  resistance  to  internal  thermal  0.1 

resistance,  I<"R” 


the  system  to  maximize  output  power  and  therefore  overall 
system  efficiency  defined  in  Eq.  (4);  Fig.  3  shows  the  profile  of 
the  current  in  accordance  with  Eq.  (10c)  for  the  thermoelectric 
module.  Fig.  4  provides  the  corresponding  temperature  profiles 
along  the  module.  These  figures  indicate  that  the  current  is  great¬ 
est  where  the  largest  temperature  differences  exist  across  the  leg 
pairs.  For  this  particular  system,  a  maximum  power  of  390  W  is 
reached  with  an  overall  system  efficiency  of  3.4%  as  defined  by 
Eq.  (4).  This  overall  efficiency  is  significantly  less  than  the  leg  pair 
efficiency  of  6.7%,  defined  by  Eq.  (5).  A  similar  simulation  for  the 
co-current  flow  configuration  generates  very  similar  results  owing 
to  the  high  heat  capacity  rate  of  the  coolant,  (mcp)c  . 

A  reduction  of  the  heat  capacity  ratio,  (mcp)c/(mcp)h,  from  ten  to 
three  and  one,  results  in  smaller  temperature  differences  and 
therefore  power  generation,  given  by  350  W  and  261  W,  respec¬ 
tively.  In  the  latter  case,  the  temperature  difference  across  the  legs 
and  current  is  nearly  constant  across  the  entire  system. 

Fig.  5  shows  that  the  addition  of  thermoelectric  modules,  each 
consisting  of  127  leg  pairs,  to  a  system  will  not  always  increase 
power.  In  fact,  for  the  parameters  of  Tables  2  and  3,  the  power 
peaks  at  around  490  W  for  approximately  65  modules.  Beyond  this 
amount,  an  increase  in  the  number  of  modules  results  in  a  negligi¬ 
ble  increase  in  generated  power  due  to  a  decrease  in  the  average 
temperature  difference  in  the  system  leg  pairs  over  the  length  of 


Fig.  3.  Leg  current,  I,  along  the  length  of  system  in  accordance  with  Eq.  (10c). 
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Fig.  4.  Temperature  profile  for  the  current  profile  shown  in  figure.  (A)  Th,  (B)  Ti,  (C) 
T2,  and  (D)  Tc. 


the  system.  For  the  situation  where  the  heat  capacity  rate  of  the 
coolant  is  three  times  that  of  the  exhaust  heat  capacity  rate,  the 
output  power  decreases  with  the  addition  of  modules  at  some 
point.  Presumably,  this  occurs  because  adding  more  modules 
results  in  parasitic  heat  losses  due  to  conduction  through  the 
added  modules.  The  peak  power  for  this  reduced  coolant  heat 
capacity  rate  is  400  W  with  40  modules. 

The  relative  resistances  of  the  system  heat  sinks  have  a  signifi¬ 
cant  impact  on  overall  system  performance.  The  results  presented 
above  were  for  the  case  where  the  ratio  of  thermal  resistances  for 
the  heat  sinks  on  both  the  hot  and  cold  sides  of  the  system  to  the 
module  thermal  resistance,  ICR^  and  I<"R''}  are  0.1.  Increasing  the 
resistance  ratio  to  0.5  while  holding  all  other  variables  constant 
(again,  those  in  Tables  2  and  3)  results  in  154W  generated,  a 
60%  decrease  in  power  generated.  This  result  indicates  that  close 
attention  is  needed  to  the  design  of  the  heat  sinks  in  any  waste 
heat  recovery  systems. 

Although  ZT,  the  thermoelectric  figure  of  merit,  is  often  used  to 
compare  potential  thermoelectric  materials,  it  is  not  sufficient  in 
open  heat  recovery  systems  consisting  of  many  modules.  For 
example,  if  the  electrical  resistivity  is  decreased  and  the  thermal 
conductivity  is  increased  by  a  factor  of  two  and  all  other  system 
parameters  are  held  fixed,  the  output  power  decreases  from 
390  W  to  381  W  despite  the  fact  that  ZT  is  unchanged.  Similarly 


Table  3 

Ranges  of  system  configurations  considered. 

(mcp)cl(mcp)h 

#  Modules  in  system 
I<"R" 

kp  &  kn  (W/m  K) 


1,  3,  10 
1,  5,  25,  50 
0.01,  0.1,  0.5 

1.1  &  1.3,  1/3*(1.1)  &  1/3*(1.3) 


if  we  increase  the  electrical  resistance  and  decrease  thermal  con¬ 
ductivity  by  a  factor  of  two,  the  output  power  generated  decreases 
to  320  W  even  though  the  ZT  is  the  same  for  both  cases.  Therefore, 
ZT  alone  is  not  a  sufficient  metric  for  comparing  materials  options 
for  open  heat  recovery  systems. 

4.1.  Alternative  optimization  approaches 

Three  alternative  approaches  have  been  used  in  the  literature  to 
optimize  systems  having  multiple  thermoelectric  leg  pairs,  and  it  is 
judicious  to  compare  results  with  those  of  the  theoretical  limit 
model.  The  three  approaches  used  are:  (1)  to  optimize  local  effi¬ 
ciency  across  each  thermoelectric  leg  pair  [19],  (2)  to  fix  the  volt¬ 
age  across  all  leg  pairs  to  optimize  the  overall  system  power  by 
adjusting  the  single  voltage  across  each  leg  pair  [18],  and  (3)  to 
string  all  legs  in  series  and  optimize  power  for  a  single  operating 
current  [16,17].  The  latter  two  approaches  are  simpler  to  imple¬ 
ment  than  the  former  in  real  systems.  None  of  these  alternatives 
will  necessarily  guarantee  the  maximum  system  efficiency  as 
defined  by  Eq.  (4).  What  follows  is  the  formulation  of  the  equations 
required  for  the  system  optimization  using  these  modeling 
approaches,  followed  by  a  comparison  with  the  theoretical  limit 
model  developed  in  Sections  2  and  3.  All  notation  used  is  as 
defined  previously. 


4.2.  Optimization  by  local  efficiency 


Following  the  approach  of  Refs.  [5,19],  the  efficiency  across  a 
single  leg  can  be  optimized  by  the  following  expression  for  current 
flux 


a(Tt  -  T2) 
R"(l+m) 


(20a) 


where m  =  v/T+zT;  z  =  and  T  =  Tl± III 
K  Ka  2 


(20b) 


We  augment  this  cited  work,  however,  by  including  thermal  resis¬ 
tance  between  the  fluids  and  thermoelectric  junctions,  which  have 
already  been  included  in  Eq.  (11)  and  illustrated  schematically  in 
Fig.  2.  After  substitution  of  jmax  into  Eq.  (11),  the  results  are 


=  gtT\ 


(Ti  -  T2) 
R"(l  +  m) 


+  K”(JA 


t2) 


«2(Ti  -r2)2 
2R"(1  +  m)2 


(21a) 


(h  -  T2) 
R"(l  +  m) 


+  I<"(h 


T2)4*ZizlJ 

2K"(l+m)2 


(21b) 


Eq.  (21)  are  divided  by  K"  and  rewritten  in  terms  of  z  according  to 
Eq.  (20b)  to  yield  the  following 


Tft-T, 

D2 

t2-tc 

Ds 


=  zTi 


=  zl2 


(Ti  -  T2) 
(1  +  m) 
(Ti  -  T2) 
(1  +m) 


+  (T,  -  T2)  - 


+  (T,  -  T2)  + 


z(T i  -  T2)2 
2(1  +m)2 
zjTi  -  T2)2 
2(1+ m)2 


(22a) 

(22b) 


Fig.  5.  Optimal  power  versus  number  of  modules  in  system  for  three  different 
heat  capacity  ratios  where  (A)  (mcp)c/(mcp)h  =  10,  (B)  (mcp)c/(mcp)h  =  3,  and  (C) 

(mcp)c/(mcp)h  =  1. 


where 

D2  =  RpC  and  D5  =  R"K" 


(22c) 
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The  set  of  nonlinear  Eq.  (22),  coupled  with  the  differential  Eqs.  (13a) 
and  (13b)  and  their  corresponding  boundary  conditions  for  co-  or 
counter-current  flows  (given  respectively  by  Eqs.  (19a)  and  (19b)), 
are  well  posed  to  solve  for  the  temperature  field  and  the  current 
flux.  We  solve  the  system  using  a  full  Newton’s  iteration  in  a  similar 
way  to  that  described  in  Section  3.  Once  solved,  the  power  is  deter¬ 
mined  as: 

p  =  w  /o  (q'^  -  q")dx  =  wf  ~Tl  -  dx  (23) 

The  model  described  above  will  be  referred  to  as  the  Local  Efficiency 
Optimization  Model. 


4.3.  Optimization  for  fixed  voltage 


For  the  case  where  the  voltage,  V ,  is  the  same  for  each  leg  pair, 
the  current  flux  is  expressed  as 


J  = 


oc(h-T2)  V 


K  K 

The  total  power  can  then  be  extracted  as 

wV2L 


w  [ LjVdx  =  ^  [\h-T2)dx- 

Jo  Ke  JO 

maximum  power  i 


K 


The  maximum  power  is  achieved  when  =  0,  which  yields 


(24) 


(25) 


(26) 


An  expression  for  current  flux,  jmax,  that  maximizes  the  power  is 
obtained  after  substitution  of  Eq.  (26)  into  Eq.  (24)  to  yield 


jmax  =  ¥[(h-T2)-r} 


(27a) 


where  V  =  — 


hi V.-w* 


(27b) 


Upon  substitution  of  Eq.  (27)  into  Eq.  (11),  dividing  by  K",  and  sim¬ 
plifying,  the  result  is 


Th  -  h 

D2 

T2-Tc 

d5 


(T1-T2)+z(r  +  ?Pj[(h-T2)-r}  (28a) 

(r1-r2)+z^T-^[(r1-r2)-r]  (28b) 


(V2  _ 

where  z  =  —77—77 ,  D2  =  R^K",  D5  =  R"I< "  and  T  given  in  eq.  (20b) 
K  Re 

(28c) 


Eq.  (28),  coupled  with  the  differential  Eqs.  (13a)  and  (13b)  and  their 
corresponding  boundary  conditions  for  co-  or  counter-current  flows 
(given  by  Eqs.  (19a)  and  (19b),  respectively),  are  well  posed  to  solve 
for  the  temperature  field  and  the  current  flux.  Once  solved,  the 
power  is  extracted  by  combining  Eqs.  (25)  and  (26)  to  yield 


P  = 


wLV2 

~K~ 


(29) 


where  V  is  given  by  Eq.  (26).  The  system  is  solved  using  a  full 
Newton’s  iteration  in  a  similar  way  to  that  described  in  Section  3. 
For  purposes  of  later  discussion,  we  denote  the  above  model  as 
the  Constant  Voltage  Optimization  Model. 


4.4.  Optimization  via  constant  current 


used  to  optimize  the  power  with  respect  to  current  flux.  This  is 
obtained  by  imposing  a  current,  solving  the  system  (13),  calculat¬ 
ing  the  power  (14),  and  then  repeating  this  process  by  systemati¬ 
cally  changing  the  current  to  find  the  highest  power.  This  model 
is  referred  to  as  the  Constant  Current  Optimization  model  in  what 
follows. 

4.5.  Comparison  of  optimization  approaches 

The  temperature  and  current  profiles  for  the  parameter  values 
in  Tables  1  and  2  are  shown  in  Figs.  6  and  7.  The  profile  for  the  the¬ 
oretical  limit  model  developed  in  Sections  2  and  3  is  also  included. 

For  this  particular  case,  the  system  efficiency  for  the  theoretical 
limit  model  is  only  0.1%  greater  than  that  of  the  local  efficiency 
optimization  model,  while  the  local  efficiency  optimization  model 
has  a  0.3%  greater  module  efficiency  as  defined  by  Eq.  (5).  This  neg¬ 
ligible  difference  was  common  for  the  wide  range  of  system  config¬ 
urations  defined  in  Table  3. 

Differences  increase  between  the  local  efficiency  optimization 
model  and  the  theoretical  limit  model  as  the  ratio  of  hot  side 
external  resistance  to  internal  thermal  resistance,  K"R thermal 
conductivity,  /<;  and  number  of  modules,  N,  decrease.  For  a  case 
where  K!'RJ'h  =  K"R!'C  =  0.01  and  k  is  decreased  by  a  factor  of  three 
from  the  base  case,  the  local  optimization  model  yields  a  total  sys¬ 
tem  efficiency  that  is  2.7%  below  that  of  the  theoretical  limit 
model.  With  the  exception  of  systems  with  one  and  five  modules, 
K"R h'  =  K"R!'C  =  0.01  and  k  is  decreased  by  a  factor  of  three  from  the 
base  case,  the  differences  between  the  two  models  were  less  than 
3%  and  in  most  cases  less  than  1%.  No  case  defined  by  the  permu¬ 
tations  in  Table  3  resulted  in  differences  greater  than  8%. 

The  differences  between  the  theoretical  limit  model  and  the 
voltage  optimization  model  are  more  pronounced  and  can  be  quite 
significant  in  some  cases.  For  the  parameter  values  in  Tables  1  and 
2,  the  constant  voltage  optimization  model  results  yields  nearly  6% 
less  power  than  that  from  the  theoretical  limit  model.  This  differ¬ 
ence  increases  dramatically  as  more  modules  are  added  to  the 
system  as  shown  in  Fig.  8  or  when  operating  the  system  in  a 
co-current  configuration.  For  example,  a  system  with  50  modules 
will  deliver  16%  less  energy  for  a  constant  voltage  condition  versus 
the  theoretical  limit  for  a  counter  flow  configuration;  by  contrast,  a 
co-current  system  with  50  modules  for  the  same  conditions  results 
in  a  nearly  a  28%  difference.  Differences  are  also  more  significant 
for  lower  mass  flow  rates  and  when  a  large  portion  of  the  thermal 
energy  is  removed  from  the  exhaust  stream.  Although  it  is  simple 


x/L 


For  a  fixed  constant  current,  Eq.  (13)  are  solved  directly,  and  the 
power  is  extracted  from  Eq.  (14).  An  inverse  parabolic  method  is 


Fig.  6.  Current  profile  along  the  length  of  system.  (A)  theoretical  limit  model,  (B) 
local  efficiency  optimization  model,  (C)  constant  voltage  optimization  model,  and 
(D)  constant  current  optimization  model. 
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Fig.  7.  Temperature  profile  for  the  current  profile  shown  in  Fig.  5.  (A)  Th,  (B)  Tlt  (C) 
T2,  and  (D)  Tc.  Solid  lines  are  for  theoretical  limit  model  and  dotted  lines  are  for 
constant  voltage  optimization  model.  Local  efficiency  optimization  model  results 
are  indistinguishable  from  theoretical  limit  curves. 


#  of  modules 

Fig.  8.  Power  versus  number  of  modules  for  counter  flow  systems  with  constant 
voltage  and  constant  current  models  compared  to  the  theoretical  limit  with  system 
parameters  defined  in  Tables  1  and  2.  (A)  theoretical  limit  and  local  efficiency 
optimization  approach  (indistinguishable  on  scale  of  plot),  (B)  constant  voltage 
approach  and  (C)  constant  current  approach. 


to  implement  a  single  fixed  voltage  across  the  entire  system  in 
practice,  the  loss  in  overall  performance  can  be  significant. 

Fig.  8  shows  that  optimal  power  predictions  can  differ  signifi¬ 
cantly  between  the  theoretical  limit  model  and  the  constant  cur¬ 
rent  optimization  model  as  the  number  of  modules  increase; 
however,  the  differences  are  much  more  pronounced  for  the  con¬ 
stant  voltage  optimization  model.  For  the  parameter  values  in 
Tables  1  and  2,  the  constant  current  optimization  model  results 
yields  more  than  2%  less  power  that  from  the  theoretical  limit 
model. 

5.  Conclusions 

The  theoretical  limit  of  the  maximum  power  that  can  be 
extracted  from  a  thermoelectric  system  in  contact  with  a  hot 
exhaust  stream  has  been  determined  using  a  variational  method. 
The  limit  is  determined  from  a  computationally  efficient  optimiza¬ 
tion  algorithm  for  co-  and  counter-current  flow  configurations 
used  in  waste  heat  recovery  systems.  The  theoretical  limit  model 
predicts  that  there  is  an  optimum  number  of  thermoelectric 


modules  that  maximize  the  power  extracted  for  any  system,  and 
that  adding  more  modules  beyond  this  optimum  can  actually 
degrade  system  performance.  Three  different  electric  loading 
approaches  found  in  the  literature  were  compared  to  the  theoret¬ 
ical  limit  using  typical  parameter  values  corresponding  to  automo¬ 
bile  exhaust— optimization  of  local  thermoelectric  efficiency, 
optimization  of  power  for  fixed  current,  and  optimization  of  power 
for  fixed  voltage.  It  is  found  that  in  most  cases,  optimization  of 
local  thermoelectric  efficiency  leads  to  predictions  close  to  that 
of  the  theoretical  limit;  deviations  arise,  however,  in  situations  of 
high  heat  transfer  resistances  in  the  fluids.  For  cases  where  little 
thermal  energy  is  extracted  from  the  exhaust  gases,  power  optimi¬ 
zation  approaches  employing  either  constant  voltage  or  constant 
current  assumptions  predict  optimal  system  power  near  the  theo¬ 
retical  limit.  For  cases  where  a  significant  portion  of  the  thermal 
energy  is  utilized  to  generate  power,  wiring  thermoelectric  leg 
pairs  in  series  (constant  current)  results  in  higher  system  efficien¬ 
cies  compared  to  cases  where  the  leg  pairs  are  wired  in  parallel 
(constant  voltage).  This  work  demonstrates  that  the  thermoelectric 
figure  of  merit  alone  is  not  sufficient  to  fully  characterize  the  effi¬ 
ciency  of  a  system.  Other  system  parameters,  included  in  the  the¬ 
oretical  limit  algorithm,  are  necessary  to  fully  characterize  the 
theoretical  system  limit. 

Due  to  its  low  computational  load,  the  theoretical  limit  model 
could  be  imbedded  in  the  design  of  a  thermoelectric  system;  in 
practice  however,  its  implementation  would  require  the  use  of  a 
relatively  sophisticated  control  scheme  to  assure  predicted  current 
loading.  The  utility  of  the  model,  however,  is  not  limited  to  this 
use,  as  it  provides  a  theoretical  limit  against  which  the  power 
drawn  from  simpler  configurations  and  control  schemes  may  be 
compared.  This  allows  engineers  to  make  judgments  on  the  bene¬ 
fits  and  costs  of  simpler,  application-specific  configurations,  trad¬ 
ing  off  complexity  (and  thus  cost  and  robustness)  with  the  loss 
of  power.  For  example,  based  on  the  case  study  considered  and 
shown  in  Fig.  8,  one  can  conclude  that  the  power  obtained  from 
a  system  of  up  to  30  modules,  run  under  a  constant  current  sce¬ 
nario,  is  approximately  that  of  the  theoretical  maximum.  Thus, 
an  engineer  would  no  doubt  choose  the  simpler  constant  current 
implementation  approach  in  an  application  where  30  modules  pro¬ 
duce  the  desired  level  of  power. 
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